﻿using System;
using SoftwareConsulting.BI11.Studio.InternalServices;
using SoftwareConsulting.BI11.NumericalMethods.DhbFunctionEvaluation;
using SoftwareConsulting.BI11.NumericalMethods.Curves;
using SoftwareConsulting.BI11.NumericalMethods.Statistics;
using SoftwareConsulting.BI11.NumericalMethods.Iterations;

namespace SoftwareConsulting.BI11.NumericalMethods.Statistics
{    
    /// <summary>
    /// Класс бета-распределения
    /// </summary>
    public sealed class BetaDistribution : ProbabilityDensityFunction
    {        
        private double _shape1;//первый параметр формы распределения
        private double _shape2;//второй параметр формы распределения        
        private double _norm;//норма распределения (кэшируется для ускорения)        
        private GammaDistribution _gamma1;//гамма-распределение для shape1, чтобы генерировать случайные числа (кэшируется для ускорения)        
        private GammaDistribution _gamma2;//гамма-распределение для shape2, чтобы генерировать случайные числа (кэшируется для ускорения)        
        private IncompleteBetaFunction _incompleteBetaFunction;//незавершенная бета-функция для распределения (кэшируется для ускорения)
        
        /// <summary>
        /// Конструктор класса
        /// </summary>
        /// <param name="shape1">первый параметр формы распределения</param>
        /// <param name="shape2">второй параметр формы распределения</param>
        public BetaDistribution(double shape1, double shape2)
        {
            if (shape1 <= 0)
                throw new ArgumentOutOfRangeException("shape1");
            if (shape2 <= 0)
                throw new ArgumentOutOfRangeException("shape2");            
            SetParameters(shape1, shape2);
        }

        /// <summary>
        /// Конструктор класса
        /// </summary>
        /// <param name="histogram">Гистограмма, из которой вытягиваются параметры по наилучшим предположениям</param>
        public BetaDistribution(Histogram histogram)
        {
            if (histogram.Minimum < 0 || histogram.Maximum > 1)
                throw new ArgumentOutOfRangeException();
            double average = histogram.Average;
            double variance = histogram.Variance;
            double a = ((1 - average) / variance - 1) * average;
            if (a <= 0)
                throw new ArgumentOutOfRangeException();
            double b = (1 / average - 1) * a;
            if (b <= 0)
                throw new ArgumentOutOfRangeException();
            SetParameters(a, b);
        }

        /// <summary>
        /// Возвращает среднее распределения
        /// </summary>
        public override double Average
        {
            get 
            { 
                return _shape1 / (_shape1 + _shape2); 
            }
            set 
            {
                Assert.Fail(Assert.Id.Assert00276);//для бета-распределения недоступно             
            }
        }

        /// <summary>
        /// Задает новые параметры
        /// </summary>
        /// <param name="shape1"></param>
        /// <param name="shape2"></param>
        private void SetParameters(double shape1, double shape2)
        {
            _shape1 = shape1;
            _shape2 = shape2;
            _norm = GammaFunction.LogBeta(_shape1, _shape2);
            _gamma1 = null;
            _gamma2 = null;
            _incompleteBetaFunction = null;
        }

        /// <summary>
        /// Задает генератор случайных чисел
        /// </summary>
        private void InitializeRandomGenerator()
        {
            _gamma1 = new GammaDistribution(_shape1, 1.0);
            _gamma2 = new GammaDistribution(_shape2, 1.0);
        }
        
        /// <summary>
        /// Возвращает значение функции плотности распределения
        /// </summary>
        /// <param name="x"></param>
        /// <returns>вероятность, нахождения случайной переменной меньше или равной, чем X</returns>
        public override double GetDistributionValue(double x)
        {
            return IncompleteBetaFunction.GetValue(x);
        }

        /// <summary>
        /// Возвращает незавершенную бета-функцию
        /// </summary>
        private IncompleteBetaFunction IncompleteBetaFunction
        {
            get
            {
                if (_incompleteBetaFunction == null)
                    _incompleteBetaFunction = new IncompleteBetaFunction(_shape1, _shape2);
                return _incompleteBetaFunction;
            }
        }

        /// <summary>
        /// Возвращает эксцесс распределения
        /// </summary>
        public override double Kurtosis
        {
            get
            {
                double s = _shape1 + _shape2;
                return 3 * (_shape1 + _shape2 + 1) * (2 * s * s +
                                        (_shape1 + _shape2 - 6) * _shape1 * _shape2)
                                / ((_shape1 + _shape2 + 2) *
                                        (_shape1 + _shape2 + 3) * _shape1 * _shape2)
                                    - 3;
            }
        }

        /// <summary>
        /// Возвращает название распределения
        /// </summary>
        public override string Name
        {
            get 
            { 
                return "Beta distribution"; 
            }
        }

        /// <summary>
        /// Возвращает/устанавливает параметры распределения
        /// </summary>
        public override double[] Parameters
        {
            get 
            { 
                return new double[] { _shape1, _shape2 }; 
            }
            set 
            { 
                SetParameters(value[0], value[1]); 
            }
        }

        /// <summary>
        /// Возвращает случайное число в соответствие с распределением
        /// </summary>
        /// <returns></returns>
        public override double GetRandomValue()
        {
                if (_gamma1 == null)
                    InitializeRandomGenerator();
                double y1 = _gamma1.GetRandomValue();
                return y1 / (y1 + _gamma2.GetRandomValue());
        }

        /// <summary>
        /// Возвращает/устанавливает первый параметр формы распределения
        /// </summary>
        public double Shape1
        {
            get
            {
                return _shape1;
            }
            set 
            { 
                SetParameters(value, _shape2); 
            }
        }

        /// <summary>
        /// Возвращает/устанавливает второй параметр формы распределения
        /// </summary>
        public double Shape2
        {
            get
            {
                return _shape2;
            }
            set 
            { 
                SetParameters(_shape1, value); 
            }
        }

        /// <summary>
        /// Возвращает коэффициент ассиметрии распределения
        /// </summary>
        public override double Skewness
        {
            get
            {
                return 2 * Math.Sqrt(_shape1 + _shape2 + 1) * (_shape2 - _shape1)
                                / (Math.Sqrt(_shape1 * _shape2)
                                                        * (_shape1 + _shape2 + 2));
            }
        }

        /// <summary>
        /// Возвращает как строку
        /// </summary>
        /// <returns></returns>
        public override string ToString()
        {
            return string.Format(Name + "({0:0.00000},{1:0.00000})", _shape1, _shape2);
        }

        /// <summary>
        /// Вычисляет значение, основываясь на функции плотности
        /// </summary>
        /// <param name="x"></param>
        /// <returns></returns>
        public override double GetValue(double x)
        {
            return Math.Exp(Math.Log(x) * (_shape1 - 1) + Math.Log(1 - x) * (_shape2 - 1) - _norm);
        }

        /// <summary>
        /// Возвращает дисперсию распределения
        /// </summary>
        public override double Variance
        {
            get
            {
                double s = _shape1 + _shape2;
                return _shape1 * _shape2 / (s * s * (_shape1 + _shape2 + 1));
            }
        }
    }
}
